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Abstract 

Corrections of order a 4 Ry are calculated for the singlet states l 1 5o and 2 1 Sq of the helium atom. The 
result for l 1 5o state is in slight disagreement with that of Korobov and Yelkhovsky in [Phys. Rev. Lett. 87, 
193003 (2001)]. The results obtained lead to a significant improvement of transition frequencies between 
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I. INTRODUCTION 



In this paper we present an approach for obtaining precise energy levels of light few electron 
atoms and perform calculations for singlet states l 1 ^ and 2 1 S of the helium atom. This approach 
is based on Quantum Electrodynamic (QED) theory and relies on expansion in the fine structure 
constant a yj bfl. It allows one to systematically include nuclear recoil effects, electron self- 
energy, and vacuum polarization. Several calculations have already been performed for triplet 
states of helium, the most accurately known being the fine structure of 2 3 Pj, where all corrections 



up to order of m a 7 have recently been obtained ||3l|4|]. Other examples include energies of l 1 ^ 
iSBl, 2 3 Si 2 3 P Jsl, which have been obtained up to the order of ma 6 . For other atoms like 
lithium iQ] or beryllium [10], energy levels are less accurately known, namely up to order ma 5 , 
but are still accurate enough to allow the nuclear charge radius to be determined from isotope shift 
measurements Jill . 

According to QED theory, the expansion of energy levels in powers of a has the following 
form: 

E{a) = E® + e<$ + £ (5) + £ (6) + + 0(a 8 ), (1.1) 

where E^ n > is a contribution of order ma n and may include powers of In a. Each term E^ can 
be expressed as an expectation value of some effective Hamiltonian or in some cases of a nonlocal 
operator. E^ = E is the eigenvalue of the nonrelativistic Hamiltonian H , which for the infinite 
nuclear mass is 

Pa Z a\ \ x - a 



E^ is the expectation value of the Breit-Pauli Hamiltonian H12I1 . 
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is the leading QED contribution, which for singlet states is 
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with being the step function and 7 the Euler constant. Eq. (11.51) contains the many-electron 
Bethe logarithm In k Q defined by 

%Pa(H -E ) lnp^)] J2 b Pb, 
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The calculation of is the subject of this work. It can be represented as 



(1.6) 



E (6) = (#(«)) + 
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where if ( 6 ) is the effective Hamiltonian of order m a 6 . Its derivation is presented in the following 
section. Since individual terms in above equation are divergent we follow the approach of Ref. B5J] 
and use the technique of dimensional regularization, details of which are presented in Appendix A. 

in the above equation is therefore a Breit-Pauli Hamiltonian in <i-dimensions, the derivation 
of which is also included in Section II. For the next order term E^ we will use an approximate 
formula based on hydrogenic values LI 5]. 



II. DERIVATION OF EFFECTIVE HAMILTONIAN 



To derive the effective Hamiltonian, we follow Ref. O, consider the Dirac equation with 
electromagnetic field, and perform a nonrelativistic expansion by the use of the Foldy-Wouthuysen 
transformation, (see Appendix B for details) 
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(2.1) 



where higher order terms are neglected. This Hamiltonian defines an effective nonrelativistic QED 
theory with Lagrangian 

C = J2 ^ d t ~ H FW )*P a + C EM , (2.2) 

a 

where Cem is the Lagrangian of the electromagnetic field and the summation goes over all par- 
ticles. We consider now the equal time retarded Green function G = G({r^}, f ; {r a }, t), where 
by {r a } we denote the set of coordinates for all particles of the system. In the absence of time 
dependent perturbation G = G(t' — t). The Fourier transform of G in the time variable t' — t can 
be written as 

G < £ > - E^hwy (23) 

which is the definition of the effective Hamiltonian H cS (E). In the nonrelativistic case H cS = Hq, 
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where [Q] e is the d = 3 — 2 e extension of the operator Q at d — 3, see Appendix A for details. All 
the relativistic and QED corrections resulting from the Lagrangian C can be represented as 

1 1 1 1 „,„s 1 1 
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(2.5) 



" E — H Q — T,{E) E - H eS (E) ' 

where T,(E) is the n-particle irreducible contribution. The energy level can be interpreted as a 
pole of G(E) as a function of E. It is convenient to consider the matrix element of G between 
the nonrelativistic wave function corresponding to this energy level. There is always such a cor- 
respondence, since relativistic and QED effects are small perturbations of the system. This matrix 
element is 

(<P\G{E)\<P) = J = - J — , (2.6) 

where 



1 E - Hq - E(E) 17 E - Eq - a(E) ' 
a(E) = m(EM) + ^(0|S(E)|0„) (0„|S(E)|0) + . . . 



(2.7) 



Having cr(£'), the correction to the energy level can be expressed as 
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Since the last term in Eq. (12.81) can be neglected up to order ma 6 , one can consider only T,(E ). 
In most cases, the explicit state dependence of £ can be eliminated by the use of commutation 
relations. The only exception is the so called Bethe logarithm, which contributes only in order 
ma 5 . If we consider this term separately, the operator £ gives an effective Hamiltonian 
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H + £ = H Q + H {4) + H {5) + H {6) + ... 
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from which one can calculate corrections to energy levels. The calculation of £ follows from the 
Feynman rules for the Lagrangian in Eq. (12.21) . We will use the photon propagator in the Coulomb 
gauge: 
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and consider separately corrections due to exchange of the Coulomb Goo and the transverse Gij 
photon. The typical one photon exchange contribution between electrons a and b is: 
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where <p is an eigenstate of H and j% is the electromagnetic current for particle a. The first terms 
of the nonrelativistic expansion of j° component are obtained from Eq. (12.11) (terms involving 
coupling to A ) 
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and of the j component (terms involving coupling to A) 

/(£) = - + 7^- o*^ • (2-13) 
m 2m 

Most of the calculation is performed in the nonretardation approximation, namely one sets A; = 
in the photon propagator G^k) and j(k). The retardation corrections are considered separately. 



Within this approximation and using the symmetrization k Cl 
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which leads to 
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One recognizes that in the nonrelativistic limit G o gives the Coulomb interaction. However this 
term is already included in H , which means that this nonrelativistic Coulomb interaction has to 
be excluded from the perturbative expansion. Next order terms resulting from j° and / lead to 
the Breit Pauli Hamiltonian H BP . This includes corrections to the electric as well as magnetic 
interactions between electrons and the nucleus. Corrections to the kinetic energy and electron- 
nucleus interaction are obtained from (12.11) by setting e A = — [(Z a) /r a ] e , 
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The derivation of electron-electron interactions is as follows. The j° component gives relativistic 
corrections of the form 
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where r = r i2 = |fi — r 2 |. We have left out in the above pure Coulomb interaction between 
electrons and neglected higher order terms. The /component gives the following corrections: 
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The complete relativistic correction is a sum of Eqs. (I2.16I) . (I2.17L (12. 1 8b . 
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Since spin-orbit terms do not lead to divergences in the second order matrix element one can 
assume for them e = 0. The spin-spin tensor interaction, the 2nd term of Eq. (12.181) . vanishes 
for singlet states, as the total spin is zero, see Eq. ( 13.291 ) for the definition of a singlet state in 
ci-dimensions. Moreover, as a result of this definition one obtains a\ J a 2 — > —d(d — 1), so 
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Both Ha and if c contribute to E^ through second order contribution, 
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Below we derive the higher order term in the nonrelativistic expansion, namely the m a 6 Hamil- 
tonian, which we call here the higher order effective Hamiltonian It is expressed as a sum of 
various contributions 

= J2 SHi + H H + H R , (2.25) 



i=0,8 



which are calculated in the following on the basis of the Foldy-Wouthuysen from Eq. (I2.ll) . Similar 
derivation has already been performed for the case d = 3 in Ref. ||2j]. One can neglect here all 
spin-orbit terms and tensor spin-spin interactions, as they vanish for singlet states. 5H is the 
kinetic energy correction, the last term in Eq. d2.ll) . 
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SHi is a correction due to the 9 th and 10 th terms in H FW in Eq. d2.ll) . These terms involve 
only A , so the nonretardation approximation is strictly valid here. This correction 5 Hi includes 
the Coulomb interaction between the electron and the nucleus, and between electrons. So, if we 
denote by V the nonrelativistic interaction potential 
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and for later use, by £ a the static electric field at the position of particle a which is produced by 
the nucleus and the other particle b 
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then 5 Hi can be written as 
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SH 2 is a correction to the Coulomb interaction between electrons which comes from the 5 th 
term in H FW Eq. d2.ll) . namely 
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If interaction of both electrons is modified by this term, it can be obtained in the nonretardation 
approximation Eq. (12.151) . namely 
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where [<5Y( r )]e * s defined in Eq. (IA30I) and we use the identity which is valid for singlet states 
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<5iY 3 is the correction that comes from 7 th term in Eq. (12.11) 



{p, d t E}. 



(2.34) 



16 m 3 

To calculate this correction, we have to return to the original expression for one-photon exchange 
Eq. (12.1 II) . We assume that particle 1 interacts with electromagnetic field by this term, while 
particle 2 by nonrelativistic coupling e A and obtain 
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After performing the A; integral, and commuting (H — _E ) with e~ lk ' Tb one expresses this cor- 
rection in terms of an effective operator 
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5H4 is the relativistic correction to transverse photon exchange. The first electron is coupled to 
A by the nonrelativistic term 



— -p ■ A — a ij B ij , 
m 4 m 

and the second one by the relativistic correction, the 4 th and 6 th terms in Eq. (12.11) 
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It is sufficient to calculate it in the nonretardation approximation 
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It is convenient at this point to introduce a notation for the vector potential at the position of 
particle a which is produced by particle b 
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This correction can then be written as 

e 



E 



1 m ,J 



pi,2p a -A a + a^ a Ai 



a=l,2 

(pi + Pi) 



■Pi 



2r 1 



ri + M+^^£S 4 ^(r). (2.41) 
8 a 



One notices that in the nonretardation approximation any correction can be simply obtained by 
replacing the magnetic field A by a static field A a - We will use this fact in further calculations. 
5H 5 comes from the coupling 
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which is present in the 5 th term in Eq. (I2.ll) . The resulting correction is obtained by replacing the 
fields E and A by the static fields produced by the other electron 
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where [l/r\ = (V[l/r] e ) 2 . 
8H 6 comes from the coupling 
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which is present in 2 nd term of Eq. d2.ll ). Again, in the nonretardation approximation the A a field 
is being replaced by the static fields produced by the other electron 
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SH 7 is a retardation correction in the nonrelativistic single transverse photon exchange. To 
calculate this correction, we have to again return to the general one-photon exchange expression, 
Eq. (12.111) . and take the transverse part of the photon propagator 
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We assume that the product j\ (k) f 2 (~k) contains at most a single power of A; . This allows one to 
perform the k° integration by encircling the only pole fc° = on 3ft(fc°) > complex half plane 
and obtain 
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where k — \k\. By using the nonrelativistic form of f and performing the retardation expansion, 

(2.49) 



1 H — E _ (Hp — E ) 
k k 2 k 3 



+ ... 



Ep-Hp-k 

where the first one contributes to the Breit-Pauli Hamiltonian, the second term to E^\ and the 
third term gives 5E 7 
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This is the most complicated term in the evaluation. After /c-integration one obtains 
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where \{r % r J — 3 5^ r 2 )/r] e is defined in Eq. (IA28I) . 

5i7 8 is the retardation correction to single transverse photon exchange, where one vertex is 
nonrelativistic, Eq. (12.131) and the second comes from the 5 th term in Eq. (12.11) 
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With the help of Eq. (12.481) one obtains the following expression for 5E 8 
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In the expansion of 1/ (i?o ~ H — k) in Eq. (12.491 ) the first term vanishes because it cancels out 
with its hermitian conjugate and the second term is a correction of order m a 6 . After commuting 
(H — Eq) on the left one obtains the effective operator 5H$ 
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Hh is the high energy contribution which is given by the forward three-photon exchange scat- 
tering amplitude. It was calculated for m a 6 correction to the parapositronium binding energy in 
ifl^l . Following [6], we define d-dimensional spin wave function in analogy to this parapositron- 



ium work (see Eq. ( 13.291) ) and take the result with reversed sign 
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where £ is the Riemann ^-function, and we follow the convention that a common factor 
[(4 7r) e T(l + e)] 2 is pulled out from all matrix elements. 
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Hr is a radiative correction and its derivation requires a separate treatment. This is based on 
our former work for helium [ 17], and this result has also been obtained in Ref. |Q]. It is a sum of 
one- and two-loop contributions 
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At this point we have obtained all contributions of the order of m ft 6 . 



III. ELIMINATION OF SINGULARITIES 

The elimination of singularities will be performed in atomic units, which in ^-dimensions be- 
come little more complicated. The nonrelativistic Hamiltonian in natural units is 

Using coordinates in atomic units 

f -> (mft)- 1/(1+2e) f, (3.2) 

it can be written as 
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If one pulls out the factor m^ 2 ^'^ 2 ^ ft 2 /( 1+2e ) from H , than one will obtain the nonrelativistic 
Hamiltonian in atomic units. Similarly for the common factor in atomic units 

m (l-10 6)/(l+26) a 6/(l+2 6) ) (3 _ 4) 

is pulled out from all the terms, which corresponds to the replacement m — > 1, a — > 1. Such a 
factor will also be pulled out from H H in Eq. (12.551) . which will lead to appearance of logarithmic 
terms. 
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We will explore now the power of dimensional regularization. All integrals of the form 
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vanish identically by definition. We will use this fact in the following. Consider the matrix element 
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The integral with q 2 vanishes, so this matrix element is equal to 0, similarly 
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Let us introduce momenta P and p 
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where p a is a momentum of the electron a and consider the matrix element 
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This is equal to because the integrand <p{P, p) p is odd in p for the singlet states considered here. 
These matrix elements and the Schrodinger equation 
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Since the electron-nucleus divergences cancel out algebraically, one does not need the matrix 
element of (1/rf). However (1/r 3 ) is needed and is calculated as follows 

3\ 



C( / d d r<p\r)r 



-3+6 e 



Cf(j) z (0) / d d rr 



+ / d 3 r6 2 (r)r 



,7T 



8 d (r)) 



(3.16) 



where (1/r 3 ) = P(l/r 3 ) is defined in Eq. dl .5b . The matrix elements of 1/r 4 can be obtained 
from Eq. (I3.15I) . but it can also be calculated directly, 



C'i I tfV0 2 (r)[V(r- 1+2£ )] 2 



C 1 2 (-l + 2e) 2 2 (O) / d d rr- i+ie (l-C 2 r 1+ ' Ze j 2 + /tfV0 2 (r)r 



l+2e\2 



+ <7T() d (r)) 



1 



4 



(3.17) 



where we assume that 1/e and line + 7 are dropped. Similarly 

1 



1 / . . r l r J \ 




T 


( 5 y + 




2r \ r 2 J 




r 



+ 7T<T(r) 



(3.18) 



From Eqs. (13.161) and (13.171) one obtains the identity 



1\ _1_ / Z Z\ 1 , 

T > +Pi7T^Pi+P2^P2- [E + - + -] - + 6tt5 3 r . 
r d / 2r z lr A \ 7*1 r 2 / r^ 



(3.19) 
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Other identities which will be used are 

-nV 2 5 d (r) 



1 6' J P : P j 



-f 


( n 


-T) 
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" 1 






r 2/ 




h 2 





r 5 3 V ; 



(3.20) 
(3.21) 





























r 









(3 r* - r 2 ) 



.p* p. 







T 










f 


pl, 




















r 


e- 






r 6 2 ) 





Pl -P2 



Pl 



T 
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|p , p ,. (3rM-^r 2 ) ^ 



ri r 2 



—P ■ p — p- P 

r 

1 {P 2 Z 



+ -7Td d (r)P 2 
3 V ' 



--n5 d (r)P 2 
3 V ; 



Z Z 

n r 2 



(3.22) 



(3.23) 



Pi -P2 



r 



Z 1 

r\ r 2 r 



E 



7T 



5 d (r) 





T 
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(r 


f 2 







p • p — p ■ p 

r 



_ _ 1 _ _ P 4 PVr, ^ ^ 1 

-pi x p 2 -pi x p 2 - 1 P H 1 

r 4r r V r± r 2 r 



(3.24) 
(3.25) 

(3.26) 



-2 7T(5 3 (r)P 2 , 

We are now ready to eliminate divergences from matrix elements of 5 Hi operators. In the follow- 
ing we make replacement a\ a 2 — > — d(d — 1). To show this, we consider matrix elements of spin 
operators with the helium singlet wave function in d- dimensions. All terms with a single spin 
operator vanish. The spin-spin operators of the form a[ J a 2 k (S^ k /d — r J r k /r 2 ) vanish because 
of angular integration. Only operators of the form o % { a 2 Q fail to vanish, and we use a defini- 
tion which was implicitly assumed in Ref. |5J]. Namely, we postulate the existence of the charge 
conjugation operator C, such that 

C- 1 = C T , (3.27) 

(3.28) 



and singlet states ips are defined by 



1 



(V-slQi® Ws> = -Tt[QiCQ t 2 C- 1 ]. 
Using this definition one obtains 

(HQv>i j \Tps) = -d(d-l) (VslQlVs), 



(3.29) 



(3.30) 
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for an arbitrary operator Q. 

All 5Ei corrections are now transformed as follows. The first term, 5E = (5H ) becomes 



16 
1 

16 



(3.31) 



4 [(ViF) 2 + (V 2 V) 2 ] + 8 (E - Vf - Qpj (E - V) p\ + 3 
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pi, 
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e- 





All singular operators in E can be handled by Eqs. (13.141) - (13.261) . and all the singularities are 
identified in the form of 2 (O) / e. Next the terms 8E\ — 5E 8 are transformed in a similar way 

SEk = T ^{-4[(V 1 V) 2 + (V 2 ^) 2 ]-2 

~ |2 (e + ^-^j 4^%) + 2 (j5 + ^—^j 4rrZ5 3 (r 2 ) 

-2 + — + — J An Z 5 3 (r) - p 2 2 Att Z 5 3 (n) - pjArr Z 5 3 (r 2 )\, 
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pi 


pI 
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r 
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= — V 2 (5 tf (r)-^«5jP i ^ + ^5jpV ) 



(3.32) 
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(3.33) 



5E 4 = p\{E-V)—[^ + r — 
2 r V r z 



P> 2 -\(e + - + -\^b\r) 
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(3.34) 
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ri r 2 \ r 



(d-l) 



(3.35) 
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(3.36) 



5E* 



f 




4) 


f (d-l) 


" 1 " 
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32 
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(3.37) 



The calculation of the E 7 contribution is lengthier, and we split it into 4 parts, corresponding to 

Eq. d22B 

5E 7 = 5E 7A + 5E 7B + 5E 7C + 6E 7D . (3.38) 
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Each part contains singular operators which are handled similarly to that in Eqs. (13.141) - (13.261) 

(d-1) 



SE 7A 
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(3.39) 
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(3.40) 
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S 13 
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4r 4 4 V e 
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1 ,• (5 jk r 2 - 3r^' r fc ) k 

Pi + oPl- A -Pi 



(3.41) 
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1 _ 1 1 _ 

- Pi X p 2 - Pi X p 2 - - Pi ■ P2 

8 r 8 

^3^ y 1 ^ fy^ (5^^" 7""^" 7"^ 



1 it / 



1 

r 

+ 3 



Pi -P2 

<y^ fy3 ry^ ry^ 1 



p\p> 2 . 



(3.42) 



iyO ryO ryi^ 

At this point we have completed the elimination of singularities from the effective Hamiltonian. It 
remains to consider, however, the second order matrix element 3a 

1 



Ha 



(E - H )> 

which requires subtractions of 1/e singularities. For this we use the transformation 



so that 



Q 



H'a + {-^0 — 3 , Q}, 



z z 

T\ r 2 



+ 



(d-1) 



(3.43) 

(3.44) 
(3.45) 



3 A = 3' A + 3 A , 



E' 



A 



H' 



H' 



A (Eo-H )'" A/ - 
E" A = (Q(E -H )Q)+2(H A )(Q)-({H A , Q})=X 1 +X 2 + X 3 . 



(3.46) 
(3.47) 
(3.48) 



E' A is finite in the limit e — > 0, and 



H' 
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(E -vy-ti — ip> + 
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Z r 1 



2r 



ryL ry3 



Z r 1 



- T -- Vi---- Vi 

4 rf 4 rf 



(3.49) 
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where Vf V\ is understood as a differentiation of on the right hand side as a function (omitting 
8 3 (r)). What remains is the calculation of Xj terms. The first two are simple 



X 1 
X 2 



1 

32 



4 ~t~ 4 
L r l r 2 J 



+ 



16 



Z / f i f 2 



\ 4:\r 1 r 2 J 2 r 



(3.50) 
(3.51) 



To calculate X 3 we split it again into four parts correspondingly 

Pi Pt 



i 



x 



— + ~ 



[5 d (r 1 ) + 5 d (r 2 )] 



p{ + n(d-2) 5 d {r) 
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(d-1) 
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L J e 



3A "I" X 3B + X 3C - + X 3 £> 



(3.52) 



and calculate each part separately 



X3A = \ {{p\+pt)Q) 
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(3.53) 



X 3B = Z (Z - 2) - 5\r x ) l + Z(Z-2)j 5\r 2 ) 
4 ro 4 
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3C* 
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2\ 1 
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'1 ?"*2 if I if y 



71 r3 , . / Z Z 

-<r(r) 1 

2 V ; V n r 2 
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(3.54) 

1" 

r 

(3.55) 
(3.56) 



The last term to be considered is H H in Eq. (12.55ft . One transforms it to atomic units by the 
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replacement in Eq. (13.21) and division by factor in Eq. (13.41) . As a result one obtains 

1 . . 39C(3) 32 . . 7\ 7T 



H 



H 



39C(3) 32 , . , 
4 In a ^ + — - 6 ln(2) 



m 



H 



7T 2 7T 2 V ' ' 3 



7T 



4 



**(r) 



4 In a -5 d (r). 



(3.57) 



At this point we have separated out all singularities. They always have the form of 2 (O)/e, and 
finally cancel between themselves. The sum of all terms, which is the main result of this work, is 
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+E' H + E' A + E C + E m + E R2 - ln(a) vr^(r), 
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(3.58) 



where E' H = (H' H ) from Eq. (l3~57l ^ is defined in Eq. (l3~47l £ c in Eq. d4~T2b . £ m and E R2 
in Eqs. (I2.56I) . (12.571) correspondingly. In addition to various identities in Eqs. d3.19l) - (l3.26l) . we 
used two further equations 

r 

Z Z I 

T\ r 2 r 



Z Z 
r\ r 2 

( r * r * _|_ £*j 

pi- 5 — V 2 



-2E^ - E 



(3.59) 



+ 



P\p\ 



-irZ[5 3 ( ri ) + 5 3 (r 2 )] + 2vr5 3 (r), 



(3.60) 



to simplify final expression. The logarithmic term in Eq. (13.581) agrees with that obtained in Ref. 
1 18]. The sum of "soft" operators (Eq. (13.581) without the last line) after (1 <->• 2) simplification 
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becomes 

h Q - — 2 8~ 8 8 ¥ I 24 

h — - — ^ Qi — y ^ 8 + 4 ^ 9 ^ — ^ — + Qn — E Z q 13 

2 ^ z z z z z z 

—Z Qu + Z Qi 5 — — Qie — — Q17 + — Qis + — Q19 — — Q20 + Q21 

Z 2 Z 2 Z 1 Z E Z 

Q21 + — Q22 + Q23 + — Q24 — — Q25 — j Q26 — g- Q27 — 2 Q28 

+ ^29 + ^30, (3.61) 
where Qi are defined in Table I. 

IV. NUMERICAL CALCULATIONS OF MATRIX ELEMENTS 

The helium wave function is expanded in a basis set of exponential functions in the form of 
Hi 

0(ri,r 2 ,r) = £ ^[e^n-to^r + ^ ^ ^ (41) 
i=i 

where aj, /3j and 73 are generated randomly with conditions: 

A 1 < a, < A 2 , A + Ti > e, ( 4 - 2 ) 
fii < A < 5 2 , Oj + ji > e, (4.3) 
Ci < 7i < C 2 , a* + A > e. (4.4) 

In order to obtain a highly precise wave function following Korobov Jl^l . we use double set of 
the form (I4.ll ). Each parameters A+, B^C^e are determined by the energy minimization, with 
the condition that e > 0, which is necessary for the normalizability of the wave function. The 
linear coefficients in Eq. (14.11) form a vector v, which is a solution of the generalized eigenvalue 
problem 

Hv = ENv, (4.5) 

where if is a matrix of the Hamiltonian in this basis set, iV is a normalization (overlap) matrix, 
and E an eigenvalue, the energy of the state corresponding to v. For the solution of the eigenvalue 
problem with J\f = 100, 300, 600, 900, 1200, 1500 we use LU decomposition in quad and sextuple 
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precision. As a result we obtain the following nonrelativistic energies in au 

£0(1%) = -2.903 724 377 034 119 592(6), 
E (2 1 S ) = -2.145 974 046 054 417 311(50). 



(4.6) 
(4.7) 



These values agree with the more accurate result of Korobov [19] and of Drake in fl20j. The 
calculation of matrix elements of nonrelativistic Hamiltonian can be performed with the use of 
one formula: 



16 7T 2 



d 3 ri / d 3 r 2 



-ari — (3r2— "fr 



7"i T 2 T 



1 



(a + /3)(/3 + 7 )( 7 + a) 



(4.8) 



The result with any additional powers of in the numerator can be obtained by differentiation 
with respect to the corresponding parameter a, [3 or 7. The matrix elements of relativistic cor- 
rections involve inverse powers of ri, r 2 , r. These, can be obtained by integration with respect to 
corresponding parameter. This leads to the appearance of logarithmic and dilogarithmic functions, 
for example 



16 7T 2 
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16 vr 2 
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d 3 r 2 



d 3 r 2 



-ari —f3r2—jr 
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(4.9) 



r 7T 2 
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+Li 2 [ I 



6 2 V/? + 7 



a + 7 



Li 2 1 



a + 7 
^+7 



. (4.10) 



a + (3 j 

All matrix elements involved in the m a 6 correction, see Table I, can be expressed in terms of 
rational, logarithmic and dilogarithmic functions, as above. The high quality of the wave function 
allows us to obtain precise matrix elements of operators and the numerical results are pre- 



sented in Table I. Some of these matrix elements have already been calculated in 112 ill and results 
in Table I are in agreement with them. 

The calculation of second order corrections E' A and E c is more complicated. The spin algebra 
in the second order matrix element Ec is simplified with the help of 

|^o) ( l S \ = \S ) (S \ (l - yY (4.11) 

where (r 1( f 2 \S ) is the wave function without the spin, and s = (<Ji + <? 2 )/2, so that 

1 



E 



c 



So 



c 



C 



Sc 



A z fn _ r 2 \ 1 r , , 
° = T W >< Pi - ^3 x + ^ ^ x ( Pl +p 2 ). 



(4.12) 
(4.13) 
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TABLE I: Expectation values of operators entering H^> for IS state, r = r\ — r 2 



Qi 


= 


47r5 3 (r 1 ) 






22.750 526 


Q2 


= 


4vr5 3 (r) 






1.336 375 


Qs 
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47r5 3 (ri)/r 2 






33.440 565 


Q± 
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4n5 3 (ri)pl 






49.160 046 


Qs 
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An 5 3 (r)/ri 






5.019 714 


Qe 
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47rS 3 (r)P 2 






18.859 765 


Qi 
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l/r 






0.945 818 


Q 8 
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1/r 2 






1.464 771 


Q9 
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l/r 3 






0.989 274 


Q10 
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l/r 4 






-3.336 383 


Qn 
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l/r? 






6.017 409 


Q12 
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l/(ri r 2 ) 






2.708 655 


Ql3 
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l/(n r) 






1.920 944 


Ql4 
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l/(ri r 2 r) 






4.167 175 


Q15 


= 


V(r?r 2 ) 






9.172 094 


Qw 
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l/(r 2 r) 






8.003 454 


Qn 


= 


i/(nr 2 ) 






3.788 791 


Qis 
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(r\ • r)/(rf r 3 ) 






3.270 472 


Q19 
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(f*i • r)/(r 3 r 2 ) 






1.827 027 


Q20 
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r « r J ^ r ? r j _ 3 gij r 2-j j ^3 


r 2 r) 




0.784 425 
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14.111 960 
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= 


Pi A 2 Pi 






21.833 598 
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Pi /r 2 pi 






4.571 652 


Q24 
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p*( r i r j + r 2 )/(n r 3 )^ 




0.811 933 


Q25 
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pi ( 3r i r i _gij r 2y r 5pj 






-3.765 488 


Q26 






'//* — foij j f 
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-0.266 894 


Q27 




2 2 

P1P2 






7.133 710 


Q28 




Pi Aip! 






37.010 642 


Q29 




Pi x P2 A pi x p 2 






4.004 703 


Q30 




P1P2 rV fc /r 3 - <5 ifc 


r j r l/ r 3 + 3 r i 




-1.591 864 
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The inversion of the operator E — H in this expression is performed in a basis set of even parity 
functions with I = 1 of the form 

0(r l5 r 2 ,r) = J^fi x f 2 [ e - a ^-^~^ r + ( n <-> r 2 )] , (4.14) 

i 

The values of parameters A4, Bi and Cj corresponding to <p are obtained by minimization of E c , 
and results of these calculations are presented in Table II. The calculation of E' A is similar to that 

TABLE II: Contributions to E^ for IS and 2S states of the helium atom. Elg is the logarithmic correction, 



last term in Eq. (13.58 



m a 6 


HeC^S) 


He(2 1 S) 


AE 


Eq 


15.465 431 


12.310132 


-3.155 299 


E' H 


-0.278 403 


-0.022 641 


0.255 762 


E' A 


-18.495 345(50) 


-16.280186(10) 


2.215159(50) 


E c 


-0.392 621 


-0.033 790 


0.358 831 


Subtotal 


-3.700 937(50) 


-4.026 485(10) 


-0.325 547(50) 


Eri 


141.924 288 


100.971873 


-40.952 415 


Er2 


1.144012 


0.890 559 


-0.253 453 


Elg 


1.643 823 


0.133 682 


-1.510141 


Total 


141.011185(50) 


97.969 630(10) 


-43.041555(50) 


-E D (Re+) 


4.000 000 


4.000 000 


0.000 000 


-E m (Re+) 


-97.971914 


-97.971914 


0.000 000 


-E m {Ee+) 


-0.873 699 


-0.873 699 


0.000000 


E^(Re) -£( 6 )(He+) 


46.165 572(50) 


3.124017(10) 


-43.041555(50) 



of E c , but additionally requires a subtraction of the reference state from the implicit sum over 
states. We obtain this by orthogonalization of H' A \(f)) with respect to eigenstate with closest to 
eigenvalue of H — E. This eigenvalue is not exactly equal to 0, because we use a basis set 
with different parameters, which are obtained by minimization of E' A . Results are presented in 
Table II. Surprisingly, the total nonlogarithmic exchange contribution after subtraction of He + 
value is very small, namely 0.299063 for IS state and —0.026485 for 2S state. This contribution is 
much smaller than the dominating one-loop contribution E R1 which is 43.952 374 and 2.999 959 
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correspondingly. It is similar for triplet states 2 3 S'i and 2 3 P\ of helium, and means that higher 
order corrections can be well approximated by the one-loop self-energy contribution. 



V. SUMMARY 

We have derived the complete order m a 6 contribution to energy levels of singlet states of 
helium. It is expressed as the expectation value of the operators in Eq. (13.581) . A similar, but not 
identical set of operators have been obtained previously by Yelkhovsky, Eq. (97) in [6], and the 
results obtained here are in slight disagreement (see Appendix C for details). The matrix elements 
of operators entering Eq. (13.581) for the ground state of helium atom are presented in Table I, 
and a few of them are in disagreement with the results presented in Ref. J5J]. Because of these 
discrepancies, calculations presented here should be verified before definite conclusions can be 
made. 

In this work we performed numerical calculations for the ground l 1 ^ and excited 2 1 S , states, 
and the results are presented in Table II. While the calculation of the Eq operators was quite 
complicated, their contribution to E^ is relatively small. The dominating contribution comes 
from the one-loop electron self-energy Em, and is given by Dirac delta functions, see Eq. 92.561 . 

The summary of all known contributions to l 1 ^ — 2 1 S transition is presented in Table III. 
The nonrelativistic energy here, E^ is the sum of fiE from Eqs. (14.61 ) and (14.71 ) with \i being the 
reduced mass, and mass polarization corrections from Ref. 



taken from Ref. 
nuclear recoil [126 



Relativistic contribution E^ is 



2511 and includes nuclear recoil. Leading QED contribution E^ also includes 
, and we use Bethe logarithms as obtained by Drake in Ref. n2n\ . Our value 
for E^ is greater by about 12 MHz from that of Korobov and Yelkhovsky in Ref. ||5j], which 
is 40 483.98(5) MHz, and we do not understand a reason of this discrepancy. We have not been 
able to find in the literature the result of Drake for E^ as well as separate results for higher order 
terms. However, the total result of Drake, see Table III, is in agreement with that of Ref. [5]. 

is obtained here and our result is greater by about 25 MHz from the result in Ref. which 
is 834.9(2). The source of this deviation is explained in Appendix C. E^ includes all electron- 
nucleus terms of order m a which are known from the hydrogen Lamb shift (one-, two- 
and three-loop contributions) and are extended to helium in the standard way. Our value is larger 
by 12 MHz from the result of Ref. [5] —84(42), because we include all a 7 terms, not only the 
leading In 2 (a). The current theoretical uncertainty comes mainly from the approximate treatment 
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TABLE III: Contributions to IS and 2S ionization energies of the helium atom in MHz. Physical constants 
from [22], Ry= 10973731.568525(73) m" 1 , a = 1/137.03599911(46), /X e = 386.1592678(26) fm, 
m a /m e = 7294.2995363(32), r a = 1.673 fm, c = 299792458. The uncertainty for E^ is due to its 
approximate calculation and is roughly estimated to be about half of E^ 





u(l'S) 


v{2 1 S) 


^u{2 1 S-l 1 S) 




-5 945 262 288.61 


-960322 874.90 


4984939413.71 


£(4) 


16 800.32 


-11974.80 


-28 775.12 


£(5) 


40495.81 


2 755.14 


-37 740.68 


£(6) 


861.24 


58.28 


-802.96 


E (7) 


-72.(36) 


-4.(2) 


68.(34) 


E-ps 


29.58 


1.99 


-27.59 


Theory 


-5 945 204174.(36) 


-960332 038.(2) 


4984 872136.(34) 


V.K. and A.Y. [5] 


-5 945 204 223.(42) 






Drake [20] 


-5 945 204 223.(91) 


-960 332 041.(25) 


4984872182.(91) 


Exp. [24J 


-5 945 204 356.(48) 


-960332 041.01(15) 


4984 872 315.(48) 


Exp. [22] 


-5 945 204238.(45) 







of these higher order terms. The exact calculation of is at present very difficult, due to high 
complexity in the derivation of and thus limits the accuracy of theoretical predictions. Our 
final theoretical predictions are in moderate agreement with the measurement of Eikema et al [23], 
and disagree significantly with the measurement by Bergeson et al in [24]. 

Having the exact formula for m a 6 corrections for singlet and as well as triplet states [Q], it is 
possible now to improve theoretical predictions for higher excited states of helium, as well as light 
helium-like ions. The extension of this approach to three- and more electron atoms or molecules 
is possible, but not all technical problems in calculating matrix elements in explicitly correlated 
basis set have been resolved yet. 
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APPENDIX A: DIMENSIONALITY REGULARIZED QED OF BOUND STATES 

The dimension of space is assumed to be d = 3 — 2 e. The surface area of <i-dimensional unit 
sphere can be obtained by considering the following integral 



1 = J d d ke~ k \ (Al) 
In Cartesian coordinates it is a product of d one-dimensional integrals 

- /> -id 

1= dke~ k2 = n d/2 , (A2) 

while in spherical coordinates it is 

I = J dn d dk k d ~ x e~ k2 = ft d X - r(d/2). (A3) 

From comparison of Eq. (IA2I) and (IA3I) one obtains 

2 7r d / 2 

* = W2V (A4) 

The d-dimensional Laplacian is V 2 = d % d\ For spherically symmetric functions /, g 

d d rV 2 {f)g = - J d d rV(f)-V{g) = -n d J drr d - l d r fd r g 

= n d J drd r (r d - l d r f)g = J d d r r 1 " 1 d r {r d ' 1 d r f) g , (A5) 

Laplacian takes the form 

V 2 = r 1-d d r r^dr. (A6) 

The photon propagator, and thus Coulomb interaction preserves its form in the momentum repre- 
sentation, while in the coordinate representation it is 

vw = /(0S e " >= ^ 1/2r(1/2 - £) '' 2< ""^ <A7) 

The alternative derivation of V(r) which omits the calculation of Fourier transform is the follow- 
ing. Consider the equation 

V 2 V(r) = -Ait8 d {r). (A8) 
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If one assumes that V(r) is of the form V(r) = C r 7 , then for r ^ 

r^drr^driCr') = 0, 



(A9) 



and therefore 7 = 2 — d = 2e — 1. The coefficient C is obtained by considering the integral with 
the trial function / 

4tt/(0) 





r d d rV 2 {V)f = 


J d d rV{V)-V{f) 


= lim 


J drr d - l tt d d r (y)d T (f) 


= lim 




— I drQ d d r (r 

r=e J £ 


= lim VL d E d ~ l d £ (C e 2 - 




= n d (2-d)Cf(o), 





therefore 



C = Ci 



An 



7T e " 1/2 r(l/2- e ). 



(A10) 



(All) 



(d-2) n d 

We are now ready to consider quantum mechanics in (^-dimensions The nonrelativistic Hamilto- 
nian of the hydrogen-like systems is 

■ 2 d 



-Z 



.1-2 £ 



(A12) 



and of helium-like systems 



H - ^ 4- ^ 7 Cl 7 Cl + Cl 



(A13) 



1 '2 '12 

The solution of stationary Schrodinger equation H = E we denote by 0, and will never refer 
to its explicit (and unknown) form in rf-dimensions. Instead, we will use only the generalized 
cusp condition to eliminate various singularities from matrix elements with relativistic operators. 
Namely, we expect, that for small r = r 12 



0(r) ^0(O)(1-CV 7 ), 



(A14) 



with some coefficient C and 7 to be obtained from the two-electron Schrodinger equation around 

r = 

[-V 2 + V(r)]0(O) (1 - CV 7 ) ^£0(O)(1-CV 7 ). (A15) 
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From the cancellation of small r singularities of the left hand side of above equation, one obtains 



7 

C = C 2 



l + 2e 
1 



4 



7r £ - 1/2 r(-l/2-e). 



Therefore, the two-electron wave function around r 12 = behaves as 

0(f 1; f 2 ) « 0(r 12 = 0) (1 -C 2 ri+ 2e ). 



(A16) 
(A17) 

(A18) 



Apart from the Coulomb potential V(r) in the coordinate space, we need also other functions, 
which appear in the calculations of relativistic operators, namely 

d d k 4tt 



V 2 (r) 
V 3 (r) 



(2 7l) d fc 4 

d d k 4tt 



(2 7r) rf k 6 

They can be obtained from the differential equations 



- V 2 V 2 M 



V(r), 



-V 2 V 3 (r) = V 2 (r), 



with the results 



V 2 (r) = C 2 r 1+2e , 
V 2 (r) = C 3 r 3+2e , 



with C 2 defined in Eq. (IA17I) and 



C3 = ^vr- 1 / 2 r(-3/2- e ). 



(A19) 
(A20) 

(A21) 
(A22) 
(A23) 

(A24) 
(A25) 
(A26) 

(A27) 



Using Vi we calculate various integrals involving photon propagator in the Coulomb gauge, 
namely 

<l"k 47T _ , ( , ,-r = 5V y 2 + 



(2vr) d A; 4 



7r e-l/2 r -l+2e 



k 2 

^ ^ T(-l/2 - e) r 2 + ~ r(l/2 - e) r* r J ' 



— (rV -3<f J r 2 ) 
8 r 



(A28) 
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and 



and 



d d k 4tt / ,, tik? 

5 J 

(2ir) d k 2 V k 2 



i k-r 



5 ij V + d i d j V 2 



^ r(l/2 - e) r 2 + T(3/2 - e) r< r 3 



2r 3 



r 2 + r* r J ) 



d d fc , / ... k { ti 

- — 4 7T 5 y — 

(2vr) d V fc 



(A29) 



a \ a 



^_il^4 7T^(r) +7r -V2 r -5+2e 

a 

- <f j 4 77 5 3 (r) + -1 (3 r l r j - 5 ij r 2 ) 



-2 5 ij T(3/2 - e) r 2 + 4 T(5/2 - e) r l r j 
= St 



(A30) 



APPENDIX B: FOLDY-WOUTHUYSEN TRANSFORMATION IN d-DIMENSIONS 



The Foldy-Wouthuysen (FW) transformation (2J, |29J is the nonrelativistic expansion of the 
Dirac Hamiltonian in an external electromagnetic field. Here we extend this transformation to the 
arbitrary dimension d of space. The Dirac Hamiltonian in the external electromagnetic field is 



H = a-7T + pm + eA 



where n = p — e A, 





( a* 1 




fl \ 


a' = \ 






10-/] 




io* J 





and 

{a\a j } = 2 6 ij I. 

The FW transformation S 12911 leads to a new Hamiltonian 



# FW = e iS {H-id t )e 



-i s 



(Bl) 



(B2) 



(B3) 



(B4) 



which decouples the upper and lower components of the Dirac wave function up to a specified 
order in the 1/m expansion. Here we calculate FW Hamiltonian up to terms which contribute to 
m a 6 correction to the energy. We use a convenient form of Foldy-Wouthuysen operator S, which 
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can be written as 



S = \ P a ■ k - ^ p (a ■ k) 3 + ^{a ■ k , eA° - id t ] + -^(a ■ 

2m l 6m z 2 m 5 m 4 



5 



a ■ E H [a • ff, [a • n , a • E]] { (a ■ if) 2 , a ■ E\ > . (B5) 

2Am A 3m* 1 J J 



4 m 2 

The F.W. Hamiltonian is expanded in a power series in S 



6 



ftw = ^H (j) + ... (B6) 

3=0 



where 



= H, 

n {1) = {iS,h {0) -id t ], 

H ij) = -[iS, n (j - 1] } for j > 1, (B7) 

and higher order terms in this expansion, denoted by dots, are neglected. The calculations of 
subsequent commutators is rather tedious but the result simple 



Hfw = eA° + ^^-^—^ + ^—^ r - — [a-7r,a-E]---- 3 {7r,d t E} 
2 m 8m 6 lorn 8m z 16 m- 1 



[a-7fJa-7rJa-7f,a-E]]] + -^ 1 {(a-7f) 2 , [a • vf , a ■ £] ) . (B8) 

lorn 4 I J 



128 m 4 

There is some arbitrariness in the operator S, what means that H FW is not unique. The standard 
approach I29L which relies on subsequent use of FW-transformations differs from this one in 
d = 3, by the transformation S with some additional even operator. 

Our aim here is to obtain the Hamiltonian for further calculations of m a 6 contribution to energy 
levels of an arbitrary light atom. For this one can neglect the vector potential A in all the terms 
having m 4 and m 5 in the denominator. Moreover, less obviously, one can neglect the term with 
a ■ A a ■ E and the B 2 term. It is because they are of second order in electromagnetic fields 
which additionally contain derivatives, and thus contribute only at higher orders. After these 
simplifications, i^FW takes the form 

Hfw = eA° + ^-^a ij B^-^--^(v-E + a ij {E\7r j }) 
2 m Am sm 6 om z \ J 

+— e —{a i *B i i , p 2 } e — {p, d t E\ + -^— A {a ij E'pP , p 2 } 

' lP 2 ^ A °]]-7^rAp 2 ^ 2A °} + ^ CB9) 



128m 4 n 64m 4 ' J 16m 
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where 



° ij = o% (BIO) 

B ij = d l A j -d j A\ (Bll) 
E i = -V i A°-d t A i . (Bll) 



APPENDIX C: COMPARISON TO FORMER WORKS 

A similar, but not identical set of operators in Eq. (13.581) has been obtained by Yelkhovsky, Eq. 
(97) in |Q]. When his operators are transformed by using the following three equations 

^(-■P-l)) = ](pi\pi+P2\p2-^(e + - + -) +4n5 3 (r)\, (CI) 
r 6 \r 2/ / 4 \ r z r z r z \ n r-i J / 

U^.p l + Z]\ = i/^lp 1 -l^ + ----^V^^ 3 (r 1 )V (C2) 



PiP\ ( -—rz- ~ + 3 ) M 



r\\r\ J I 2 \ rf r i V T\ r 



g*W( -3 — + 3 ^^)PlP2 



~ P* 3rVt P j + t W - ^ - 5 3 (r) P 2 , (C3) 

8 r° 4 12 

then almost all operators agree with one exception. The difference between operators in Eq. (97) 
of Ref. U and that of ours Eq. (l3~58l) is 

A = ^n5 3 (r)P 2 . (C4) 

For this reason we checked the calculation in Ref. [6]. The derivation of initial operators was very 
similar to our former work in [ 17]. However, the electron-electron Coulomb interaction, according 



to Ref. [5] involves the term 

Ip-Api+Pim, (cs) 

while our calculations in fl 1 VII give 

The difference between Eq. (IC5I) and Eq. (IC6I) for singlet S-states is equal to A from Eq. (IC4I) . 
and this term should be subtracted from Eq. (97) of Ref. Jg]. Although in this work we use a 
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different formalism we obtain the same result (after using Schrodinger equation) as in Ref. 0. 
namely the sum of our terms Vq N + V^ e + 5qV differs from 5 Hi + 5H 2 + 5H 3 by exactly the 
same term A. 

pi 

Considering numerical matrix elements, we found that results presented in Ref. yD are not reli- 
able. Most of them are accurate to only three digits, for example the matrix element corresponding 
to Q28, according to KY is 36.983, while our result is 37.010 642. Some of them are accurate only 
to the first digit, for example for Q 29 KY gives 4 x 1.078 = 4.312 and our result is 4.004 703. Some 
of matrix elements in Ref. Jfl contain misprints in their presentation, namely they should include 
additional 1/2 on the left hand side, to agree with numerical values and to be consistent with Eq. 

n 

(97) of Ref. 1 6]. Most importantly, some matrix elements from [5] are in error, for example 



2r 3 \r r 2 

By using integration by parts we transform this matrix element to the form 

3 % ( r l\ \ 3/1 1 \ 

-•P-o ) = 7W-^- 47r <K r )> = -4.246525, (C8) 



2r 3 \r 2J / 4 \ r 4 r 3 
which is in disagreement with the result of J5|], namely —0.958. Alternatively, we use the identity 
in Eq. (13.191 ) and obtain the same numerical value as in Eq. dC8l) . Another example is for 
which result of [5] is 2 x 4.749 and our result from Table I is —0.266 894. In conclusion, the 
numerical matrix elements of Ref. Jf] should be verified. 
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